Open in Colab

Cognestic 2023 Functional Connectivity in Python


Petar Raykov & Rik Henson

This tutorial is using the nilearn and nltools toolboxes to install them run the cell below

Why Python:


Comprehensive standard library

Exceptional external libraries

For some other useful courses on neuroimaging in Python see:

  1. General overview
  2. DartBrains full undegraduate course
  3. Naturalistic fMRI data analysis course
  4. Brainiak tutorials, a toolbox for advanced fMRI analyses
  5. CBU Video tutorial

Functional Connectivity (FC) in fMRI

Often we are interested in examining how certain brain regions are activated by a task. However, a different way to analyse fMRI data is to examine how different regions may be intrinsically connected. Regions can be structurally connected or functionally connected. Here we will go through how to perform functional connectivity analyses

The simpliest way to do this is to perform Pearson correlation of the time-series of one region - A with the time-series of another region - B. This should be interpreted as temporal covariation between regions as this does not necessarily imply that the two regions are directly communicating with each other, see figure below and lecture.

To compute direct covariation one can perform partial correlation, which takes into account the connectivity between two regions whilst regressing out the common effects from other regions. For example, it shows the connectivity between A and B after regressing out C from both A and B.

Image

The most common way to compute functional connectivity however remains using pearson correlation. To perform this correlation analysis however, first we have to pre-process the data and extract the regions' time-series.

In this tutorial we will show how to perform:

  1. Seed to whole brain Connectivity:
  2. Perform ROI to ROI connectivity
  3. Do a Mass-Univariate approach to examine which edges correlate with Age
  4. Perform a Classification analysis using whole ROIxROI matrix
  5. Extract time-series from regions and compute single-subject ICA
  6. Perform Group ICA

First we load some data

nilearn comes with datasets that have been made freely available and can be easily downloaded. Alternatively data can be downloaded through openneuro fetch_openneuro_dataset or with datalad - see neuroimaging example.

We load a single subject first. The datafetcher is described here and the data is cited here. This is an FMRI data while children and adults were watching a short movie. However, we would treat it as if participants were doing a resting-state scan for now.

Here we load the confounds table and the image path and image

Here we load the Functional Image for a single subject

Lets for now drop the aCompcor columns and compute 24 motion parameters from the 6 motion paramers

Let is create a function to clean up the confound files

Simple image manipulaitons

You can skip this section as it has been covered previously, however if you want a refresher you can have a look at the code below, we calculate the global signal that can be used later on as an additional confound




A useful way to visualise the data is with a carpet plot introduced in 1

Note that this visualization is very helpful for understanding NiftiMaskers and how to represent 4-D fMRI data as matrix on which we can easily perform mathematical operations


see also 2 for some applications of the plot. 3 for a recent discussion on global signal regression

The NiftiMasker

It transforms the data from 4D to a 2D matrix (time by voxels) this makes it easy to work with various functions

This can be used to compute brain mask but also can be used to temporally pre-process your data and extract data from ROIs.

Importantly one can transform the 2D matrix back into 4D image








1. We first show how to perform seed to whole brain analysis

Before continuing with the seed to whole brain analysis we show that one can comput also the seed to seed connectivity matrix

Now we return to the multi subject data set to show how to perform the seed to whole brain connectivity analysis

  1. First we have to extract the seed time course from each subject
  2. We have to extract the voxel time courses from the masked fMRI image from each subject
  3. Then we can compute the correlation between the seed and every voxel in the brain
  4. Then we save the seed to whole brain Nifti images that can be submitted to second level analysis similar to what you learned before.
  5. We show how to compare the connectivity between the seed and the whole brain between children and adults.

Design Matrix or X is what we use in GLMs to define the regressors or tests we are interested in doing. You can see that the above design matrix corresponds to a two sample t-test. Note there are multiple ways one can define a design matrix to perform effectively the same hypothesis test and they affect the interpretation of the beta values and what contrast is needed. see this videos for more detailed examples, specifically Day 10 video.

Mass Univariate

see here for some worked out examples of different designs

To define a contrast we can either use the column name in the design matrix or provide a numpy array

Note Nilearn does not perform the typical cluster level correction. Below we show voxel-wise FWE correction which is very strict.

In nilearn you can do cluster level inference following the method described here

Although we cannot perform parametric cluster correction that relies on random field theory. We can perform multiple correction using permutation testing. This can achieve cluster size correction or perform the TFCE method that aims to combine information from both the height and extend of the cluster when determining whether a cluster is significant.

The non_parametric_inference returns the negative log10 of the p_values. So above we are plotting negative log p-values but not thresholding them to be below alpha of 0.05, to do this we have to set the threshold at -np.log10(0.05)

Alternatively you can use fsl's randomize


fslmerge -t PCC_4D.nii $(imglob pcc*gz)
Glm # set design and contrasts with Gui
randomise -i PCC_4D.nii.gz -o FSLTwo_sample -d fsl_design.mat -t fsl_design.con -m EPI_Mask.nii -T -n 5000
fsl_design.mat is the design matrix like dummy_coded fsl_design.con is a list of contrasts e.g., for child and child - adult. T is for TFCE n is number of permutations

For theory of cluster correction you can see this video 1 and about TFCE this one

2. We show how to extract signal from an Atlas and compute ROIxROI connectivity Matrix

To choose ROIs we may use an Atlas or Parcellation

Parcellations can be anatomically, or functionally defined.

Two types of atlases can be fetched from nilearn automatically:

  1. labeled images 3D with each ROI being different number; No overlap across ROIs
  2. probibalistic 4D image which each ROI in the 4th dimension and there is overlap

Here we show a labeled or 3D nifti image that holds each region as a different number. There is no overlap across ROIs

Here we load a probibalistic atlas that holds each ROI as a separate volume in a 4D image. ROIs can overlap. This is similar to how we used the ICA components as ROIs

There are multiple atlases and possible parcellations. This is still an ongoing research topic and it is unlikely that there is one parcellation that fits all questions. Below I include some important references on the topic. It is likely that there is not a single parcellation that works best for all cases. Your choice of an atlas will be guided on what you want to do later on

Check out Neuroparc and the paper Lawrence et al., (2021)

You can create anatomical or functional atlases see Moghimi et al., (2021) for a recent review, see also Passingham et al., (2002)

The numer of ROIs to be used is an ongoing research question using higher number of ROIs seems to better capture the underlying functional connectivity present at the voxel scale see this. The functional parcellation could also differ based on the task has been previously shown Geerligs et al., (2015) see also this and this and possibly this; and this

Extracting signal from probalistic atlas - NiftiMapsMasker

The nifti masker here extracted the average ROI signal, but we may want to standartise the average time courses, furthermore we may want to apply a temporal filter and regress out certain confounds. These things can be achieved with the NiftiMasker

We use Connecitivy Measure module

it returns a N x K x K array representing N -subjects - here 1; K x K connectivity; K is number of regions

**Note**

given the 2D array you can also compute the correlations manually just using numpy however, note that the Connectivity measure uses not a maximum likelihood estimation of the covariance that most people are familiar with but a Ledoit-Wolf estimator that basically shrinks all correlations towards 0 see here; Below we show how to compute correlations manually both with numpy and scipy since this may be more flexible

Regress Confounds before computing connectivity

We have quite a lot of regions begin correlated with each other

However we did not regress out motion confounds here we examine how removing the confounds affects the connectivity matrix

What if we also regress out global signal

Global signal regression introduces quite a bit of negative connectivities and is still a debated pre-processing option in the field see here. Lets plot the functional connectivity matrix on a brain






Lets compute the ROIxROI over multiple subjects

Now we have for each subject their ROIxROI connectivity stored in the 40 x 741 array (39*38)/2 - We ignore the diagonal as this is self-connections.

The mass univariate approach fits a GLM model to each edge (cell in the ROIxROI array) independently. We then can examine which edges are effected by our experimental manipulation or relate to some trait of the participants.

Here we correlate each edge with Age to try and examine which edges show different connectivity across age.

Mass Univariate

We can apply an FDR correction across all the ROI connections

Alternatively we can apply FWE correction with permutation

To get at more direct connectivity we can use partial correlation rather than pearson correlation. To achieve that we can use the inverse Sparce Covariance Matrix in order to get the direct connections. Note this is similiar to computing the partial correlation it can allows us to view the association between two regions after account it for other regions influence.

For detailed example of using the Sparce inverse covariance see here

There are different Connectivity measures, full correlation, partial correlation - connectivity between to regions conditioned on connectivity of other regions

Apart from running mass univariate GLM on the connectivity edges we can also run multivariate methods which take as regressors the functional edges and try to predict some characteristic of the participants. For instance we might use the whole connectivity matrix to examine whether it containts information that can help us classify children vs adults.

We run a classification of participants age to see which correlation method better. For different classification examples check sklearn documentation here.

Note for imbalanced classes read this paper

















Here we show how to run ICA on an individual subject's data

Nilearn does not have a specific module for ICA and just uses a the sklearn ICA module that can be applied to any numpy array data see paper that introduces nilearn

If you are classifying manually components as noise check out this paper

Note for invidual subject data you can run ICA AROMA easily from fMRIprep. Other common approaches are GIFT matlab toolbox and Melodic FSL toolbox.

We have the data as a 2D numpy array. Data is treated as n_samples x n_features

Now lets use nltools which provides a better interface for ICA

we need to load the data as BrainData object see here

A problem with single-subject ICA is that the components are not necessarily consistent across participants.

One can concatenate the time-series and run ICA on the whole group of subjects this approach concatenation approach is common and used in both GIFT and FSL

If the time-course is similar across subjects e.g., they performed the same task one can use tensorial ICA as used in FSL


To keep the components consistent across subjects one can also perform Canonical ICA described in detail here This analysis takes a 3 step approach of:

  1. First reducing the dimension of each participants fMRI data;
  2. Performing Canonical Correlation analysis, a correlation measure between two sets of variables, that ensures consistency across participants and finds a lower dimensional space common to the group of participants
  3. Performing ICA on the common group space to extract group ICA components

Lets use the function we wrote to clean up the confounds from fMRIprep

see here for discussion on which confounds to use, see also here, Opinions are far and wide on this but often you will need to regress more than the 6 motion parameters and may consider including signal from white matter, CSF and possibly global signal. Including their derivatives and squares is an additional common strategy. Note ideally you will regress all confounds in a single model

Alternatively one can run Dictonary learning see below

Another common way to perform Group ICA is FSL's dual regression approach see here and lecture here

ICA or clustering approaches can be useful for identifying networks and testing hypothesis at the network level, however, one might be particularly interested how one specific regions' connectivity with the rest of the brain differs across participants. This functional connectivity analysis is often termed seed to whole brain analysis. We will below show you how to extract signal from a single seed and compute seed to whole brain connectivity

Exercises

  1. Add Spike Regressors to the clean-up pipeline
  2. Use a separate atlas to examine differences in the functional connectome
  3. Run Analyses with more people

You can check out Rik's lecture on effective connectivity and Dynamic causal modelling. We also include a matlab tutorial on DCM.